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Abstract 

I We consider a thin film consisting of two layers of immiscible liquids on a solid horizontal (heated) 

Oh ■ 

• ' substrate. Both, the free liquid-liquid and the liquid-gas interface of such a bilayer liquid film 

Ch ; 

may be unstable due to effective molecular interactions relevant for ultrathin layers below 100 nm 
thickness, or due to temperature-gradient caused Marangoni fiows in the heated case. Using a 
^ ' long wave approximation we derive coupled evolution equations for the interface profiles for the 

o\ ■ 

■ general non-isothermal situation allowing for slip at the substrate. Linear and nonlinear analyses 
O , 

, of the short- and long-time film evolution are performed for isothermal ultrathin layers taking into 



o 



X 



account destabilizing long-range and stabilizing short-range molecular interactions. It is shown that 
the initial instability can be of a varicose, zigzag or mixed type. However, in the nonlinear stage 



^ . of the evolution the mode type and therefore the pattern morphology can change via switching 

> 



between two different branches of stationary solutions or via coarsening along a single branch. 
PACS numbers: 68.15.+e, 81.16.Rf, 68.55.-a, 47.20.Ky 
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I. INTRODUCTION. 



Instabilities of thin liquid films between a solid substrate and a gas atmosphere have 
attracted much scientific interest. The main focus lies thereby on front instabilities of mov- 
ing contact linesii^iS or on instabilities of the free liquid-gas interface of a fiat filn4^. 
A recent review can be found in Ref. 7. To analyze such instabilities a long-wave or 
lubrication approximation^'^'^ is often used as a very powerful tool especially for low 
Reynolds number film flows. At present the basic behavior of one-layer films in the phys- 
ically different thickness ranges is well understood. Several instability mechanisms ex- 
ist that by means of different driving forces may destabilize an initially flat film. They 
are described, analysed and modelled in a large number of experimenta l^i^i^'^i^^i^^i^^i^^ and 
theoretica l7ii^i^ii7ii«iiQi2n,2i,22,23,24,25,26,27,28,29,sn,si ^Q^j^g^ thicknesses d less than 

about 100 nm, effective molecular interactions between the film surface and the substrate 
dominate all the other forces, like thermo- and soluto-capillarity or gravity, and thus de- 
termine the film stability. For heated films of thicknesses above 100 nm, eventually ther- 
mocapillary forces become the most important influence leading to an instability caused by 
large-scale Marangoni convection^>2S. It is dominant up to an upper limit of the film thick- 
ness determined by the competition between large-scale and small-scale convection modes^. 
For even thicker films with thicknesses above 100 fim also the gravity force becomes impor- 
tant. Depending on its direction it may stabilize the large-scale Marangoni instability or 
destabilize the film further (Rayleigh- Taylor instability) The lubrication approximation 
is valid up to a limiting film thickness obtained by the requirement that the wave length 
of the dominant instability mode Am is much larger than the film thickness d, i.e. Xm ^ d. 
For the Rayleigh- Taylor instability Am depends on the interfacial tension, the density of the 
liquid and the gravitational acceleration but not on the film thickness''. It is of the order of 
10^ to 10^ /im implying an upper limit for the film thickness of 10^ to 10'^ fim. 

The stability and evolution of liquid films of thicknesses below 100 nm is determined 
by effectiv molecular interactions between substrate and film arising, for instance, from 
Van der Waals, electrostatic or entropic interactions^!^. Such films are lineary unstable 
if the energy of the intermolecular interaction is a convex function of the film thickness. 
For film thicknesses above 10 nm the long-range Van der Waals forces dominate. They can 
be of different nature depending on the molecular properties of the involved media. One 
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distinguishes interactions between two randomly orienting dipoles (orientation interaction), 
between a randomly orienting dipole and an induced dipole (induction interaction), and 
between a fluctuating dipole and an induced dipole (dispersion interaction). Between two 
parallel interfaces at a distance rf, all these forces decay as A/(P where A is the Hamaker 
constant^. An unstable situation corresponds to a positive Hamaker constant. Note, how- 
ever that different schools use different sign conventions. The dominant wave length of the 
instability Am increases monotonically with d as A ~ rf^ (see Ref.bj). The stability of a film 
may change dramatically for a substrate coated with a layer of different dielectric properties 
as, for instance, a silicon substrate (Si) coated with an silicon oxide layer (SiO)25.. There, 
for an oxide layer of about 2nm only ultrathin polystyrene (PS) films below 4nm thickness 
are linearly unstable. Increasing the film thickness in the linearly unstable range, the wave 
length A increases rapidly and diverges at the critical thickness dc- For d > dc the film is 
linearly stable, but may rupture due to finite disturbances. 

Imagine one replaces the (solid) coating layer by a liquid layer transforming the system 
in a two-layer liquid film. Some of the results obtained for a solid coating can be directly 
transferred to the new situation. The stability of the (now) upper layer still depends on 
the (now liquid) coating layer. However, additionally the liquid coating layer itself may be 
unstable making a re-evaluation of the stability necessary. This thought experiment leads 
quite naturally to the extension of the well studied one-layer systems to two-layer systems 
that are the subject of the present work. In general, there exist two possible two-layer 
geometries. On the one hand the two liquid layers can be situated between two solid plates 
leaving only the interface between the two liquids free to move. In consequence such a system 
can be described by a single evolution equation^^'SL^. On the other hand the two layers can 
be situated between a solid substrate and a gas atmosphere. Then both, the liquid-liquid 
and the liquid-gas interface are free to move and their evolution has to be described by 
coupled evolution equations. Models were derived, for instance, assuming a lower liquid 
layer that is much thicker than the upper layep^, and for two-layer systems with surfactants 
(and non-newtonian behaviour )^2*^i*^ or including evaporation^^*^^^. A two-layer system 



46. In Ref. 



47 a similar 



under the solely influence of molecular interactions is studied in Ref. 
system is studied, however, the evolution equations are given in terms of variations of an 
energy functional. 

The experimental interest in two-layer liquid films is up to now mainly focused on the 
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dewetting of a liquid layer from a very thick layer, i.e. a liquid bulk substrat o^^i^^i^°i^^i^^ . In 
contrast Ref.lssI studies the dewetting of a polystyrene (PS) layer of 15 to 68 nm thickness 
from a 46 nm thick polyamide (PA) layer. The substrate is a silicon (Si) wafer covered with 
a layer of native oxide. At high temperature (195°C) and small thicknesses (15...35nm) the 
PS layer is unstable and dewets exhibiting typical spinodal patterns. At low temperature 
(115°C ) th e PA layer is solid resulting in a stable PS layer, independent of its thickness. 



Ref. 



5J studies relatively thick layers (100 nm to 1 fim thickness) of poly(dimethylsiloxane) 



(PDMS) layers on a liquid substrate of fluorinated PDMS. They show that the PDMS 
films are metastable and may dewet by nucleation of holes. The velocity of the growth of 
holes depends on the viscosity and thickness of the substrate. In another system (PS layer 
on poly(methylmethacrylate) (PMMA) layer, both with thicknesses of about 100 nm) the 
dewetting velocity was found to exhibit a minimum as a function of the viscosity of the lower 
layer—. Furthermore, for a polycarbonate (PC) layer on a poly(styrene-co-acrylonitrile) 
(SAN) layer Ref.0 reports a non-trivial dependence of the dewetting velocity on both layer 
thicknesses. 

Recently we presented coupled evolution equations for two-layer liquid films taking into 
account long-range Van der Waals interactions only^. This allows to study pathways to- 
wards rupture but can not describe the long-time evolution of such films as, for instance, 
necessary for the description of the above mentioned experiments. In the present paper we 
present the derivation of the system of evolution equations for a general interaction energy 
and a non-isothermal situation and then study the short- and long-time evolution of two- 
layer liquid films incorporating long-range and short-range interactions. Thereby the main 
focus lies on ultrathin layers with respective thicknesses below 100 nm for which the effective 
molecular interactions between the four media are the dominant influence. 

The paper is structured as follows. The problem is formulated in Section III Al followed 
by the derivation of the relevant long wave evolution equations for a general non-isothermal 
case allowing also for slip at the substrate. Focusing on the isothermal case without slip 
these equations are analyzed and integrated numerically in the subsequent sections. Sec- 
tion presents the linear stability analysis for fiat films starting from the general case in 
Section IIII A[ focusing on long-range van der Waals interactions in Section IIIIB[ discussing 
the possible mode types in Section IIII CI and introducing important limiting cases in Sec- 
tion I^^Dl In Section HVl we investigate the non-linear behaviour discussing in Section FlV Al 
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non-uniform stationary solutions as extrema of the Lyapunov functional, presenting in Sec- 
tion ^^B] two different pathways for mode-changes in the course of the time evolution, and 
discussing long-time stationary solutions for different experimentally relevant systems in 
Section II V CI Concluding remarks follow in Section |3 Expressions for the surface tension 
gradients in terms of gradients of the layer thicknesses are discussed in an Appendix. 



II. DERIVATION OF COUPLED LARGE- AMPLITUDE EVOLUTION EQUA- 
TIONS 



A. General case 



First, we derive coupled evolution equations for the profiles of the liquid-liquid interface 
hi{x, y) and the liquid-gas interface h2{x, y). We assume that the layers are thin enough that 
convective terms can be neglected. Considering a two-dimensional geometry as sketched in 
Fig. d the respective Stokes equations for the two layers are 

V(Pi + 0i) = /ijAtTi, (1) 

where i = 1,2 denotes the respective layer. For each layer Vi = {ui,Wi) is the velocity field, 
Pi the pressure, (pi the potential of the bulk forces and yUj the viscosity. The constant mean 
film thicknesses are denoted by di = {J^ hidx)/L and d2 = {J^ h2dx)/L where L is the 
lateral system size. A lubrication approximation^ is applied assuming the ratio of vertical 
and horizontal length scales to be small. As smallness parameter we use the ratio e = di/X 
where A is the characteristic lateral length scale of the film instability. In zeroth order in e 
the Stokes equations ((H) simplify to 

/i29^M2 = d^p2 (2) 
d,P2 = (3) 

jJid^ui = d^pi (4) 
d.pi = 0, (5) 

where the pi stand for pi + (pi. At the substrate {z = 0) we use a Navier slip and a no- 
penetration condition, i.e. 

Ml = PdzUi and Wi = 0, (6) 



respectively. The slip length is denoted by p. At the liquid- liquid interface {z = hi) we 
use the continuity of the velocity field, the kinematic condition and the continuity of the 
tangential component of the liquid stress tensor 



Ui = U2, Wi = W2, (7) 

wi = dthi + uidxhi, (8) 

(9) 

and 

/ii^^ui - /i2<9^M2 = d^ai2, (10) 

respectively. The normal stress condition is discussed below. At the liquid-gas interface 
{z = only the kinematic condition and the continuity of the tangential component of the 
liquid stress tensor apply, i.e. 

W2 = dth2 + U2dxh2, (11) 
H2dzU2 = d^a2. (12) 

The (Ji2 and a2 stand for the interfacial tensions of the liquid-liquid and of the liquid-gas 
interface, respectively. The boundary conditions for the normal component of the stress 
tensor are written incorporating the disjoining pressures at the liquid-liquid Ili{hi, /12) and 
at the liquid-gas Il2{hi, /12) interface, respectively. They represent effective molecular inter- 
actions between the interfaces that result, for instance, from Van der Waals interactions^. 
They are discussed in detail below. For the liquid-gas interface {z = /12) we obtain 

P2{h2) -po = -a2dlh2 + n2(/ii, h2) 

and for the liquid- liquid interface {z = hi) 

Piih) - p2ihi) = -ai2dlhi + Ui{hi,h2), (13) 

where po is the constant pressure in the gas atmosphere. Eqs. ()13|) can be written in terms 
of variations of an energy functional F[hi, h2] 

6F 

Pi{hi) -P2{hi) = — 
ohi 

5F 

P2{h2) - Po = jj^, (14) 



with 



CTl- 



(15) 



2 2 

and f{hi, being the free energy of the fiat films per unit area. 

Equations (j21) and (0} are integrated three times with respect to z to obtain the stream 
functions defined by {wi = —d^^^i, Ui = dz'^i). The six x-dependent integration con- 
stants are determined using the boundary conditions (jHI), 0, (fTn|l . and (fT^ . Thus the 
velocity fields in the two layers are given by 



1 1 

ui = —{d^pi) — + —{z + p)Ki 
Hi 2 /ii 

U2 = — (9^P2) — H K2{z-hi) + ui{hi), 

/i2 2 /i2 /i2 2 



(16) 



with Ki = K2 + d^ai2 + [(<9^P2) - (f^xPi)] /^i and K2 = d^a2 - {d^P2)h2- 

The stream functions \E'j are related to the flow in the lower layer Fi = J^^ ui dz and to 
the one in the upper layer r2 = /^^^ U2 dz by Fi = \l/i(/ii), Fi + r2 = ^2ih2)- Using the 
we rewrite Eqs. dSJ and (fTT|) to obtain the evolution equations for the two interface profiles 



dthi + d,[^i{hi)] = 0, 
dth2 + d^[^2{h2)] = 0. 



(17) 
(18) 



Written in terms of the energy functional they read 
dthi = d, 
dth2 = d, 

with the mobility matrices of the pressure terms 



5F 5F 

Q^^^^TT + ^i2'9^T7 Diid^ai2 - Di2d^a2 

ohi dh2 

6F 6F 

Q2ldxYl- + Q22dxYT D2ldxO-i2 - 0220^^2 

ohi dh2 



(19) 



hi 



+ I3hl 



^ih2-!f)+f3hh2 



^ [^ih2-^)+Phh2 ^^^(J^-l) + f 



(20) 



2 V ^ 3 

and of the tangential stress terms 



fii [ h,{h2 -!f)+ (3h2 + hi{h2 -'f)+ f3h2 



(21) 



2 / ' f^'"^ 2fi2 ' '"J-V"-^ 2 

respectively. Note, that the mobility matrix Q is symmetric and all mobilities Qik and Do. 
are positive. Dropping the terms representing the effective molecular interactions, Eqs. ()19|) 
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represent the fully nonlinear equivalent for the weakly nonlinear equations derived in Ref. 



55, 



56l . Assuming that the interfacial tensions are influenced by thermocapillarity only, one can 



express the derivatives dx(Ti2 and dx(T2 in terms of gradients of local thicknesses d^hi. This 
is done in the Appendix. 

For isothermal ultrathin liquid films one has {dx(yi2 = dx<^2 = 0). The situation is 
then relaxational (or variational), i.e. Eqs. ()19p possess a Lyapunov functional, namely the 
energy functional F, which decreases monotonously in time as shown next. The total time 
derivative of the Lyapunov functional is dF/dt = f (^f^dthi + j^dth2^ dx. Expressing dthi 
by Eqs. (fTI?|l and using partial integration with periodic boundary conditions, one obtains 

f^-/5«40«)- 

Because 

det Q = ih^hia + A_ht{h2 - + P (^ih2 - hr + hj^^) > (23) 
9/ii/i2 12/if \3fi( 3jUi/i2 / 

and Qii > 0, the quadratic form in Eq. ()22|) is positive definite implying dF/dt < 0. The 
existence of F allows to identify the stationary solutions of Eqs. ^] with the extrema of F. 
This will be used below in Section llV Al 



B. The disjoining pressures nj(/ii,/i2) 

In many important for instance, for polymer films on apolar substrateg^*-^^, the 

interaction energy is mainly determined by its long-range apolar dispersion part. However, 
if the model only takes into account a destabilizing long-range interaction the time evolution 
definitively leads to rupture of the upper or lower layei^ making it impossible to study the 
long-time coarsening behaviour. To be able to study the long-time evolution one has to 
include stabilizing short-range interactions2£i2^ into the model. Although these are normally 
not included for films of thicknesses above 10 nm because they do not change the stability 
of fiat films, also for such films they become important in the non-linear stage of evolution 
when the local thicknesses become comparable to their interaction length. 

The long-range part of the interaction energy for each pair of interfaces (see Fi 
resulting from dispersive Van der Waals interactions is given by Aijj^i/ 271 h"^ (see Ref. 
where Aij^i is a (four-index) Hamaker constant which corresponds to the interaction between 
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the interfaces i —j and k — l. Each index in Aijki can be one out of (yf, 1, 2 and s, denoting gas, 
hquid 1, hquid 2 and substrate, respectively. The four-index Hamaker constant is calculated 
using an equivalent of Eq. (11.13) of Ref.|33| that is based on the assumption that the main 
absorption frequencies of all involved media are about u^. = 3 x 10^^ Hz and that the zero 



frequency contribution is negligible. One uses 



{n] - n^A(n^ - nl) 



ijkl 



(24) 

where are the refractive indices of the media. The three-index Hamaker constants are 
given by 

The short-range forces which can be of an electrostatic or structural naturoi^^ decay 
exponentially with h. The electrostatic part results from the formation of diffuse electric 
double-layers in the vicinity of interfaces involving polar liquid a^^i^^i^°i^^ . For films with 
thicknesses in the range of the Debye length the diffuse double-layers at the two interfaces 
overlap resulting in a repulsive or attractive force between the interfaces. The corresponding 
interaction energy between the interfaces s — 1 and 1 — 2 is given by 5*1 exp [(/q — hi)/li] and 
between the interfaces 1 — 2 and 2 — g hj S2 exp [(/q — (/i2 ~ ^i))/^2]) where Iq = 0.158nm 
is the Born repulsion length, and h, I2 ~ L.lOnm are the interaction lengths of the short- 
range interactions^S. Further on we consider the two correlation lengths li and I2 to be 
equal and denote them by /. 5*1 > and 5*2 > are the short-range components of the total 
spreading coefficients. They are related to the lower layer on the substrate below a bulk of 
the upper liquid and to the upper layer on the lower film as substrate below the ambient 
gas, respectively. We do not take into account short-range interactions between interfaces 
s — 1 and 2 — (yf. 

Collecting the long-range and the short-range forces the disjoining pressures Hi(/ii,/;,2) 
and H2(/ii, /i2) are specified as 



Hl(/ll,/l2) 



A 



21s 



A 



12g 



Si 


k - hi 


_^S2 


— exp 

n 


[ h \ 


+ -j-exp 



h - (^2 - ^1) 



A 



12g 



+ 



A 



g21s 



6n{h2 — hiY Grrhl 



S2 

— exp 
h 



Iq - {h2 - hi) 



(26) 



To non-dimensionalize Eqs. (jl9p we scale the thicknesses with /, the lateral coordinate x 



with A = l{d2 — (ii)>/27r(Ti/|y4i2g|, and time t with r = {2n)'^ai^il{d2 — diY /A\2g- Then the 
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scaled energy functional 
F = 



+ ci{hi - di) + C2{h2 - + Si exp {-hi) + 5*2 exp {hi - h2)] dx, (27) 



involves the scaled Hamaker constants Aijki = [{d2 — di)/l]^ Aijki/\Ai2g\, spreading coeffi- 
cients Si = 2tt [{d2 — di)]^ SiBxp {do/l)/\Ai2g\, and mean layer thicknesses di = di/l. The q 
are Lagrange multipliers that ensure mass conservation for the two liquids. The correspond- 
ing energy scale is \Ai2g\/2Ti{d2 — diY and the ratios of the mean layer thicknesses, surface 
tensions and viscosities are d = ^2/(^1, cr = 0-2/0"! and fi = yU2//^i5 respectively. Further 
on, we denote the scaled variables using the same symbols as before, i.e. the scaled mean 
thicknesses are given by di and the local thicknesses by hi. The non-dimensional mobility 
matrices are obtained from Eqs. ()2()|1 and ()21|1 by dropping the factor 1/fii and replacing (3 
by ^/l. 

III. LINEAR STABILITY 

A. General stability of flat films 

We start the analysis of our model for two-layer films by discussing the linear stability 
of fiat films with hi{x) = di and h2{x) = d2- Eqs. are linearized in e ^ 1 for small 
amplitude disturbances ex^exp (7t)exp(A;x) for i = 1,2 where k, 7 and x = (x? 1) are the 
wave number, growth rate and amplitudes of the disturbance, respectively. The dispersion 
relation 7(/c) is obtained by solving the resulting eigenvalue problem {J — '~fl)x = 0. 

For the isothermal case {d^Cu = dxCr2 = in Eqs. p9|) ). the corresponding non-symmetric 
Jacobi matrix J is given by J = —k'^Q ■ E, where Q is the scaled mobility matrix. E is the 
energy matrix 

( ^ + k'^ \ 

zr I dh'j dhidh2 \ 

where f{hi, /12) is the local part of the energy density from Eq. ()27|) . This yields 



Tr /Tr^ 

7 = y±V — -Det, (29) 

where Tr = — A;^ [2Q12E12 + QuEu + (522-E22] and Det = k^ det Q det E are the trace and the 
determinant of J . Since det Q 7^ the eigenvalue problem can be written as the generalized 
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eigenvalue problem {k^E+'~fQ~^)x = 0. Because E and are both symmetric and is 
positive definite one can deduce that all eigenvalues 7 are real^ as expected for a variational 
problem. In the non-isothermal case, the Jacobi matrix is given by J = —k'^{Q ■ E — D T), 
where F is a scaled matrix of coefficients for the Marangoni terms. It is defined in Section 1X1 
Neither the matrix D nor T are symmetric. This leads to (in general) complex eigenvalues 
indicating the possibility of oscillating motion in the non-isothermal case^^^. 

Going back to the isothermal case, inspection of the generalized eigenvalue problem shows 
that the stability threshold is completely determined by the eigenvalues of E. Since the 
surface tension terms are always positive, the onset of the instability is always found for 
A; = 0, i.e. the system is linearly stable, independently of the wavelength of the disturbance, 
for 

det > and En > at k = 0. (30) 

An instability sets in if at least one of the conditions (j3U|) is violated. Then the fiat two-layer 
film is unstable to disturbances with k larger zero and smaller than a cutoff wavenumber 

'~ 2[dhl^adhl)^A\dhl adhl) ^ a\dhdh2) ^ ^ 

determined by the condition det E{kc) = 0. 

Fig.|21shows a schematic stability diagram in the plane {En, -^22)- The stability threshold 
E11E22 = -^12, Ell > is a hyperbola, represented by the solid line. The unstable region 
below and left of that line is divided by a second hyperbola into a two-mode and a one-mode 
region. In the two-mode region both growth rates given by Eq. (j29|) are positive for k smaller 
then the respective cut-off kc- In the one-mode region only one 7 is positive for k < k^. Fixing 
all other system parameters, (\.eiE{k) is determined by k. If at A; = the system is in the 
two-mode region, then by increasing k one passes two times a line det^^ = 0, as indicated 
by the dashed arrow in Fig.|21 At each crossing a growth rates (Eq. (I29p ) becomes negative, 
i.e. a mode is stabilized. If at A; = the system is in the one-mode region the line det E = Q 
is crossed only once (dot-dashed line). 

B. Long-range interaction only 

As detailed above the stabilizing short-range interaction only becomes important if at 
least one layer thickness is locally comparable to the interaction length /, i.e. if a layer 
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becomes thinner than about 10 nm. Therefore, to study the hnear stabihty of thicker layers 
one can neglect the short-range terms in Eq. (fTTj) . In this case, Eq. (jHUjl is used to study the 
role of the Hamaker constants Eq. (j24j) in the linear evolution of the system. First we note 
that the Hamaker constants are coupled through the refractive indices of the media rij. This 
allows only for selected combinations of signs of the Aij^ and Aij^i as given in Table |1J 

For fixed Hamaker constants, i.e. fixed combination of materials, det Eq = det E{k = 0) 
is a function of the ratio d of the layer thicknesses only. Using Table |l] one can show that for 
positive d'^f/dh'^ the equation det EQ^d) = can only have the solution d = 1, i.e. d2—di = 0. 
This means that only for vanishing upper layer the system can be on the stability threshold. 
In consequence the stability threshold can not be crossed by solely changing the ratio of layer 
thicknesses. This was analyzed in Ref.|4j for a variety of experimentally studied systems. 
Increasing the ratio d from one the system remains either completely in the unstable or in 
the stable region. 

To compare the stability behaviour of two-layer and one-layer films we introduce two 
effective one-layer films as follows. In (1) we assume the lower layer to be solid, i.e. we 
regard the upper layer as a one-layer film on a coated substrate. In (2) we assume the upper 
layer to be rigid but deformable by bending. The lower liquid layer corresponds then to a 
one- layer film on a solid bulk substrate. In case (1) the one-layer liquid film is unstable if the 
second derivative of the energy with respect to the film thickness /12 is negative, d"^ f /h\ < 0. 
The stability threshold at d'^f /hi = can be crossed by changing the layer thickness h2 — hi 
or the thickness of the coating layer hi. This was demonstrated in Refs. QandQ for a PS 
film on Si wafers covered with a 1.6 nm thick SiO layer. In case (2) the one-layer liquid film 
is unstable for d'^fjhl < 0. It can also be destabilized by changing the layer thicknesses, as 



was shown in Ref. 



64 for a rigid PS layer on top of a liquid PDMS layer on a Si substrate. 



Comparing the stability thresholds for the two effective one-layer systems to the stability 
diagram in Fig.|2l shows that the stability threshold of the two-layer system lies in the region 
where both effective one-layer systems are stable. This indicates that a two-layer system is 
always less stable than corresponding effective one-layer systems. 
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C. Different instability modes 



The stability threshold can be studied in rather general terms as was done above because 
its main features do not depend on surface tensions or viscosities. However, this is not 
the case for the characteristics of the instability like mode type, growth rates or dominant 
wave length. To discuss these we focus in the following on selected two-layer films stud- 
ied experiment ally^i^^^. We consider various combinations of layers of polystyrene (PS), 
poly(methylmethacrylate) (PMMA) and poly(dimethylsiloxane) (PDMS) on a silicon (Si) 
or on a silicon-oxide (SiO) substrate. The Hamaker constants for different combinations are 
calculated using Eq. (j2D) and given in Table [TTl 

The linear instability of a two-layer film has two different modes. It can be of zigzag or 
varicose type. For the former the deflections of the two interfaces are in phase whereas for 
the latter they are in anti-phase. For special parameter values one can also find a mixed 
type, where modes are present because they have equal fastest growth rates^. 

The model studied in Ref. assumes a thick lower layer thereby neglecting the interaction 
between the substrate and the liquid-liquid and the liquid-gas interface. In this case only 
the varicose mode can be unstable. In the general case, however, also the zigzag mode can 
become unstable. Both modes are normally asymmetric, i.e. the deflection amplitudes of 
the two interfaces differ. We characterize this asymmetry by = x/(l + X^)- Negative 
(positive) corresponds to varicose (zigzag) modes. The value |0| = 1/2 represents the 
symmetric case, whereas = corresponds to maximal asymetry, i.e. one of the interfaces 
is fiat. The asymmetry increases with the ratio of the surface tensions a. Note, that the 
dispersion relation and the type of the dominant mode depend on a and /i, whereas the 
stabihty diagram Fig.El does not. 

The two mode types are plotted in Fig.Elfor a Si/PMMA/PS/air system for a fixed value 
of d and different ratios of viscosities /i. The dispersion relations jik) are shown together 
with the corresponding 0. We show here that the type of the dominant mode can be changed 
by varying d or /i. It can also change in dependence of the ratio of the interfacial tensions 
a as studied in Ref.j^zl 

Strictly speaking, the concept of the mode type characterized by is only valid in the 
linear stage of the evolution. However, to discuss morphology changes we generalize this 
concept to nonlinear thickness profiles hi{x). We define a generalized mode or solution type 
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by the integral 

If {hi - di){h2 - 



-'int 



Lj [{h-diy + ih-d^)']'^''' ^^^^ 

taken over the domain length L. In many cases the sign of the product {hi — di){h2 — ^2) 
does not depend on x allowing to 'read' the mode type directly from the plots of the layer 
profiles. For small deflection amplitudes Eq. ()H2|1 gives again the linear mode type defined 
above. In the following we use the notion 'mode-type' in the linear and in the nonlinear 
regime. 



D. Limiting cases 

For general rfj the radical in Eq. ()29p does not allow to give an analytic expression for the 
wave number and the characteristic growth time = l/7m of the fastest growing mode. 
Nevertheless, one can derive asymptotic expressions for km and in the two important 
limiting cases of (1) small thickness of the upper layer d2 — di <^ d2 and (2) small thickness 
of the lower layer di -C d2. First consider case (1), which corresponds to a liquid film (the 
upper layer) on a liquid substrate (the very thick lower layer). The dimensional k^, and 
are then given by 



kr. 



\Al2g\ 



{d2- diY^ 47r(Tefr 

(33) 



16(27r)Vcfr/ii(rf2 - d-'^^ 



diA^i 



H29 

with (Jeff = 0102! [pi + 02). Note that all variables are in their dimensional form. 

Interestingly, the growth time depends only on the viscosity of the lower layer /ii 
and does not depend on [I2. This can be explained by the fact that the flow in the lower 
layer which is related to /ii, is much larger than that in the upper one^^. At constant 
thickness of the lower layer, is proportional to (^2 — c?i)^, i.e. a liquid film on a bulk liquid 
substrate evolves faster than the same film on a solid substrate (growth time ~ (^2 — rfi)^) 
and even faster than the same film on a solid substrate with slippage (growth time ~ 

{d2 - dif I\\\'i^l{d2 - di)\). 

In case (2), which corresponds to a liquid film (the lower layer) on a solid substrate below 
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the other hquid (the very thick upper layer), the dimensional km and Tm are given by 



u — 



21s 



12(27r)Vi/iirff 



(34) 



Note that in case (2) km and Tm coincide with k\Q^ and riow, respectively, the characteristics 
of the dominant mode of the instability of a liquid film below a bulk liquid calculated using 
one-layer theory. A discussion of this geometry for a Rayleigh- Taylor instability can be 



found in Refs. 



18 and 



66 



E. Long-range apolar and short-range polar interactions 

The stability analysis based only on long-range interactions becomes incorrect for layer 
thicknesses in the range of the interaction length / of short-range interactions. Practically, 
the latter become important (well) below 10 nm layer thickness. In contrast to the result for 
the exclusive action of long-range van der Waals forces, in the regime where both, short- and 
long-range interactions, are important the stability threshold can be crossed by changing the 
layer thicknesses di. Fig.E] presents a selection of qualitatively different stability diagrams 
in the plane spanned by the layer thicknesses obtained when varying the strength of the 
short-range interaction for a fixed long-range interaction. 

By changing the short-range part of the spreading coefficient 5*1 and 5*2 one finds seven 
topologically different types of such diagrams. These types correspond to regions in the 
(5*1, 5*2) plane as indicated in Fig.El In the absolute unstable region bounded on the right by 
(•S*!)!!!!!! = (e/4)^A2is/|Ai2c,| aud above by (S'2)min = (e/4)'' the system can not be stabilized 
by changing di or d2- Only if at least one of the two Si is larger than the corresponding critical 
value (a) stable region(s) can be found in the (c/i, 6/2 — ^1) plane (see Fig.|H). For Si > {Si)-[ain 
a stable region exists that extends towards infinite (^2 — c^i), as shown in Figs.lU^a), (b) and 
(d). Thereby, for large (^2 — di) the system is stable for {di)^i^ < di < ((ii)max, where 
{di)raax and ((ii)jnm are the solutions of the equation A2is/|Ai2c,| = Six'^exp (— x). Similarly, 
for 5*2 > (5'2)min & Stable region exists that extends towards infinite di, as in Figs.|H^a) to 
(d). For large di the system is stable for (^2 — c^Omin < ^2 < (c?2 — c?i)max, where (^2 — rfi)max 
and ((^2 — c?i)min are the solutions of the equation 1 = S'2x'^exp (— x). In the gray shaded 
triangle at the center of Fig. El an additional bounded stable region exists in the {di, ^2 — <^i) 
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plane (see Figs.EJ^b) and (c)). Combining the different conditions gives the following seven 
types of stability diagrams. 



I: The stable region is continuous and extends in respective stripes towards infinite di 
and — di (Fig.Efa)). 

II: There exist two separated stable regions, one extending towards infinite di and the 
other one towards infinite d2 — di (Fig.Efd)). 

Ill: Similar to Type II but with an additional bounded stable region (Fig-H^b)). 

IV: A bounded stable region exists together with an unbounded region extending towards 
infinite di (Fig.|H^c)). 

V: Similar to type IV but with the unbounded region extending towards infinite {d2 — di) 
(not shown). 

VI: Only one stable region exists extending towards infinite di (not shown). 

VII: Similar to type VI but with the unbounded region extending towards infinite (c/2 — di) 
(not shown). 

Further on we will focus our attention on the stability diagram of type I. 

IV. NON-LINEAR BEHAVIOUR 

A. Stationary solutions as extrema of the Lyapunov functional 

To find periodic stationary solutions of the scaled Eqs. ()19p. the time derivatives dthi are 
set to zero. Integration yields 



where the Ci are constants and F is given by Eq. (j2Z|). Note that the left hand sides 
of Eqs. ()35|1 represent the flow in the lower layer and the total flow, respectively. For a 
stationary state both flows are zero, i.e. the Ci = C2 = 0. Because the mobility matrix Q 




(35) 
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is non-singular, one concludes from Eqs. (jSHI) that the stationary states of the Eqs. (fT^ are 
the extrema of the Lyapunov functional F, i.e. they are solutions of 

- C^zx/il + 777- = Cl 

ohi 
df 

-adxxh2 + TTT- = C2, (36) 
dh2 

where / denotes the local part of Eq. ()27j) and the constants q correspond to the La- 
grangian multipliers introduced in Section lll Al To obtain a finite amplitude solution for 
given mean thicknesses we use continuation technique o^^i^^i^^ . We start with analytically 
known stationary periodic small-amplitude profiles, which correspond to the linear eigen- 
functions for the critical wave number kc- By continuation we follow the family of solu- 
tions changing the period L. We characterize the solutions by the deflection amplitudes 
Al and A2, the energy E, the norm L2 = {l/L) J[{hi — diY + {h2 — (^2)^] dx and the in- 
tegral mode type (pint- To determine the stability of the stationary solutions hi{x), we 
add small perturbations Shi{x) ~ exp{(3t) to both interfaces hi{x), linearize the full time- 
dependent evolution equations (fT9|l around hi{x) and solve the obtained eigenvalue problem 
L{hi, d^hi, dx)5h{x) = (36h{x) for the linear operator L after discretizing it in space. The 
sign of the largest eigenvalue {3 determines the stability of the stationary solution. Note 
that due to the translational invariance of the evolution equations (fT^ . there exists always 
a symmetry mode 6hi{x) = dxhi{x) with the eigenvalue j3 = 0. 



B. Mode type transitions 

The type of the dominant instability mode calculated above by linear stability may not 
persist in the course of the nonlinear evolution. Possible mode type changes may have a 
dramatic effect on the (observable) overall morphology of the film. We investigate these 
changes by studying both, the evolution in time of the film profiles and the stationary 
solutions obtained by continuation. 

The evolution in time is obtained by numerical simulations of the scaled coupled evolu- 
tion equations Eqs. p9p in a one-dimensional periodic domain. Both, semi-implicit pseudo- 
spectral and explicit time integration schemes are used. Initial conditions consist of flat 
layers with an imposed noise of amplitude 0.001. 
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1. Transition via branch switching 



First the time evolution of an initially flat fllm is studied for parameters as in Fig.|21(a) 
using a domain size equal to four times the fastest growing wave length A^- A time sequence 
of snap shots is shown in Fig.El In the early stage of the evolution a varicose mode develops 
{t = 8.1) as expected from the linear analysis. Then in a sub-domain of size Am the deflection 
amplitudes increase dramatically accompanied by a morphological change towards a zigzag 
type proflle (t = 9.7). This is further illustrated by the dependence of the integral mode- 
type fEq.l32j) on time given in Fig.[7|^b). Further on, the length of the zigzag part increases 
slightly, and coarsening sets in resulting in the disappearance of one varicose-type drop 
{t = 16.4). Next one of the remaining drops increases it amplitude and flips to a zigzag type 
hole (t = 20.6). Finally the last remaining varicose-type drop disappears {t ^ 28), and the 
system approaches a stationary (but not stable) state. The evolution of the relative energy 
of the proflle in time is given in Fig.[7|^a). It is seen very clearly that the phases of very 
slow evolution correspond to solutions that are close to stationary solutions. This results 
from the fact that the (unstable) stationary solutions form saddle points in function space 
that are approached along their stable manifolds and subsequently repel the system along 
their unstable manifolds (for a discussion see Ref.Q) The evolution stops after a further 
coarsening step, when the period becomes equal to the system size (not shown). 

To explain the observed mode-type change, we study the stationary solutions of the 
evolution equations Eqs. (jTH)). We flnd a family of solutions starting at the subcritical 
primary bifurcation, then turning three times at saddle-node bifurcations (folds) and going 
towards inflnite periods (see Fig.|HI^a)). A stability analysis using the solution period as the 
period of the disturbance (thereby excluding coarsening modes) shows that two branches 
are stable (solid lines) and two are unstable (dashed lines). Along the flrst unstable branch, 
which starts at and ends at the flrst fold at L ^ 60 the energy E is always larger than the 
one of the flat fllm Eq (Fig.|Hfb)), and it increases with decreasing period. This subcritical 
branch corresponds to nucleation solutions that have to be overcome to break the fllm in 
parts smaller than Lc (see Ref.0 for a discussion of this type of solutions for a one-layer 
system). The flrst stable branch starts at the flrst fold at L ^ 59 and ends at the second 
fold at L ^ 116. Its relative energy decreases monotonically with increasing period. Mostly 
it is energetically preferable to the flat fllm. The second unstable branch (between the 
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second fold at L 116 and the third fold at L ^ 79) turns back towards smaller periods. 
The second stable branch starts at the third fold and goes towards infinite periods. Its 
energy decreases rapidly from values even above the flat film to values below the ones of 
the first stable branch. The energy of the second unstable branch is always larger than the 
energies of both stable branches. This indicates that it corresponds to nucleation solutions, 
or critical solutions that have to be overcome to switch between the two stable branches. 
Along the second unstable branch, the mode- type changes from varicose to zigzag (Fig.|Hl(d)) 
explaining the non-trivial behavior observed in the time evolution shown in Fig.|Bl There 
are two stable solutions with a period equal to the dominant linear wave length (Am ~ 108) 
(see Fig. EI). The one of higher energy that is approached first in the time evolution is of 
varicose type whereas the one of lower energy that the system switches to is of zigzag type 
(cp. Fig.lHI^b)). A transition between the two solutions is accompanied by a strong increase 
of the amplitude Ai (see Fig.|Hfa)). 

2. Transition via coarsening 

A mode-type change is not always connected to a transition between different branches 
of stationary solutions. Also coarsening along one branch may lead to such a change if the 
mode-type varies along the branch. To demonstrate this, we simulate the time evolution 
using parameters as in Fig.Ofb). A time sequence of profiles and the corresponding evolution 
of the relative energy and the integral mode-type are shown in Figs. ^1 and respectively. 
Early in the evolution the layer profiles represent a zigzag mode (t = 6.0) corresponding to 
the linear results (Fig.El^b)). Then, within the very short period of time from t = 6.0 to 
10.8, nonlinear effects result in a first change towards a varicose type profile, as shown in 
the inset of Fig. llir b). Then the system approaches the branch of stationary solutions. As a 
result the evolution slows down and the pattern begins to coarse. With ongoing coarsening 
[t > 10.8) the size of the droplets increases (t = 135, t = 461) and at very late times 
(t > 490) the mode type changes back to zigzag type (Fig.HTJ^b)). Here, the amplitudes of 
the interfaces do not change dramatically, as was the case in Section IIVB II In this sense 
the transition is continuous. The characteristics of the corresponding stationary solutions 
are shown in Fig.^1 The primary bifurcation is again subcritical fFig. lT^ a)) the solution 
family continues towards smaller periods until turning at a saddle-node bifurcation (fold) 
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and heading towards infinite periods. The subcritical branch is unstable with energies higher 
than the energy of the flat film fFig. lT^ b)). The second branch starting at the fold {L ^ 26) 
consists of solutions whose energy decreases monotonically with increasing period. They are 
stable to disturbances of identical period but unstable to coarsening modes. Fig.^l(d) shows 
that the solution with the period equal to = 50.81 is of varicose type. The corresponding 
layer profile is shown in Fig. lT^ a) together with the profile after the first coarsening step. 
When the period becomes larger than 94.2 , the solution changes to zigzag type fPig. fT^ d)) 
as shown in Fig. ll3r b). This explains the mode- type change found in the time evolution 
(Fig.CH. 

Here we have restricted our attention to a parameter set corresponding to region I of 
Fig-El i-G- corresponding to the stability diagram shown in Fig.E](a). The existence of a 
stable branch of stationary solutions which continues towards infinite period implies that the 
rupture of the two layers is completely avoided by the short-range repulsion. However, this 
may not be the case for parameter ranges belonging to the other types of stability diagrams. 
A detailed analysis of the stationary solutions for all types will be done elsewhere. 

C. Large-period stationary solutions 

The stability of the numerical code, used to solve the evolution equations Eqs. (|T!I|) . re- 
quires a very small time step t = 0.00001. As a result it takes very long even to reach the 
final stationary solution in a system of size 4Am using 256 grid points. It is not feasible at 
the moment to study many coarsening steps in this way. However, one can rely on contin- 
uation techniques that use an adaptive spatial grid along the continuation patb^ to obtain 
stationary solutions of arbitrarily large periods that correspond to solutions that would be 
obtained in a time evolution at very late times. We show in Fig.El possible large-period 
stationary solutions for different physical systems that are investigated experimentally. One 
finds qualitatively different morphologies like a drop of the lower liquid 'looking through' a 
nearly fiat film of the upper layer for a Si/PMMA/PS/air system. Note however that also 
the upper layer is continuous (due to the stabilizing short-range interaction), i.e. also the 
drop is covered by a very thin layer of the upper liquid. This is more pronounced for a 
SiO/PS/PDMS/air system. In contrast, for a SiO/PMMA/PS/air system one finds a drop 
of the upper liquid 'swimming' on the lower liquid that however is attracted towards the 
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base of the drop. These equihbrium solutions are equivalents of drop configurations studied 



in Ref . 



71 



for macroscopic (but smaller than the capilary length) drops. However, here the 
mesoscopic contact angles are not given explicitely but result from the underlying effective 
molecular interactions, i.e. the short- and long-range forces used. 



V. CONCLUSION 



We have derived coupled non-linear evolution equations for the profiles of the liquid-liquid 
and liquid-gas interfaces of a thin two-layer liquid film heated from below allowing for slip 
at the substrate. We have shown that in the isothermal case the evolution equations can be 
written in terms of variations of an appropriate Lyapunov functional F which monotonically 
decreases in time. The stability conditions for flat layers have been given in terms of F. 
We have shown that a two-layer film is less stable than related effective one-layer films 
introduced in Section IIIIBl Even if both effective one-layer films are stable the two-layer 
film may be unstable if the determinant of the energy matrix det E is negative. 

We have shown that if the Hamaker constants are given by the usual expression (Eq. ()24|) ). 
i.e. they are coupled through the refractive indices, and no other forces are present, the sta- 
bility of the flat films with thicknesses of (~ 100 nm) can not be changed by solely changing 
the layer thicknesses. Incorporating a stabilizing short-range interaction the stability can 
be changed in this way. We have classified the resulting possible types of stability diagrams 
in the space of the layer thicknesses and given a 'phase diagram' in terms of the short- 
range parts of the spreading coefficients for the occurence of the different types of stability 
diagrams. 

In general, the linear stability analysis of the fiat film has shown that both, varicose or 
zigzag mode, may be unstable depending on the ratios of the layer thicknesses, viscosities 
and surface tensions (see also Ref. [4^. This seems to be in contrast to Ref.0. However, the 
difference arises because there it is assumed that the lower layer is thick compared to the 
upper layer neglecting thereby all interactions with the substrate. Then the zigzag mode is 
always stable. 

The introduction of the stabilizing short-range interaction allows to study the long-time 
evolution and stationary layer profiles. Possible stationary states have then been determined 
as extrema of the Lyapunov functional F . The resulting bifurcation diagrams show a rich 
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branch structure that depends strongly on parameter values. We have focused on one type 
of stability diagram where a stable branch going towards infinite period always exists. This 
implies the existence of a non-ruptured stationary state in the long-time limit also in the 
time evolution. 

We have found that the mode type of a profile may change during the evolution of an 
instability on three ways. First, the profile type changes in the course of the short-time 
evolution. This is connected to different mode types found for the dominant linear mode 
and the stationary solution of equal period on the solution branch approached first in the 
time evolution. It seems that this behaviour is more probable for a subcritical primary 
bifurcation. In the case studied here this change is from zigzag to varicose. In the nonlinear 
regime the profile can change its type by (i) jumping from one to another stable branch 
and by (ii) coarsening along a single stable branch. Combinations of the different ways may 
also be possible. We have found that for the parameters considered here both nonlinear 
transitions go from varicose towards zigzag type. In case (i) the transition occurs without 
change of the period, but with a dramatic increase in amplitude of the profile. In the case (ii) 
the transition occurs continuously without amplitude jump because mediated by coarsening 
a small-period varicose mode turns into a large-period zigzag one. 

In all examples considered here (except the SiO/PMMA/PS/air system with di = 30, 
d2 = 50) we have found a zigzag-type solution at large periods. For the future it would be 
very interesting to further analyze the stationary solutions for a broader range of experi- 
mentally interesting systems like the ones studied in Refs. This should clarify under 
which conditions the long-time (or large-period) solutions are energetically preferable and 
determine how 'late' the transition may occur. A systematic analysis of all types of stability 
diagrams would also discuss metastability and absolute stability of the fiat two-layer films. 
Furthermore, we are very optimistic that the evolution equations presented here will serve 
to study the questions of hole growth and possible front instabilities in the dewetting of a 
liquid layer on a liquid substrate of finite thicknes3^>^. 

APPENDIX A 

To further specify the thermocapillar part of Eqs.^Jwe rewrite the derivatives d^cru and 
dxO'2 in terms of the gradients of the local thicknesses dxhi. In the long- wave approximation^ 
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the temperature field is in both layers a linear function of the vertical coordinate i.e. 
Ti = ttiZ + hi. To determine the coefficients we consider a three-layer geometry fFig.[T3jl. 
i.e. we take into account the heat conduction in a gas layer of finite thickness dg = dt — c?2, 
where dt is the distance between the substrate and an upper plate. The temperature in the 
gas layer is Tg = agZ + bg. The boundary conditions at both interfaces are continuity of the 
temperature field and continuity of the heat flux HidzTi = HkdzTk, where is the thermal 
conductivity of the i-th layer. The temperatures at the substrate Tq and at the upper plate 
Tf are constant. The coefficients and bi depend on the local thicknesses hi and are given 
by 

aAT 



02 



ai 



where AT = Tq — T2 and 



a 



bi = To, bg = Tt- agdt 
= agKghi ( -] + To, 



d2- dt - Kghi/Ki - Kg{d2 - di)/K2 



Kghi/ni — Kg{d2 — di)/K2 
Here T2 is the temperature of the liquid-gas interface, when both interfaces are undeformed, 
i.e. for hi = di. The above formulas allow to determine the derivatives 



dxcru = TiAhi + Vi2d^h2 
dxcr2 = T2idxhi + T22dxh2, 



(Al) 



where the matrix F is determined as follows 



Kl dT ^ 



-^^h{^ - 1] 



\ 



\ ^(^c - J/.26) + Kb) - bh,{^^ - 1)} ) 



1 2 



(A2) 



Here a = (aAT) / dt — h2 + ^hi + ^{h2 — hi) , 6 = ^^(I/ki — I/K2) and c = {c/^ — /i2(l — 
^)}. For the linear normal Marangoni effect dai2/ dT and da2/dT are negative and constant. 
The Eqs. (jMl) are used in Eqs. to obtain a closed system of equations for hi and /i2. 



23 



1 A. M. Cazabat, F. Heslot, S. M. Troian, and P. Carles, Nature 346, 824 (1990). 

2 F. Brochard-Wyart and P. dc Cannes, J. Phys.: Cond. Mat 6, A9 (1993). 

3 A. L. Bertozzi, A. Miinch, X. Fanton, and A. M. Cazabat, Phys. Rev. Lett. 81, 5169 (1998). 
^ G. Reiter, Phys. Rev. Lett. 68, 75 (1992). 

^ S. J. VanHook, M. F. Schatz, W. D. McCormick, J. B. Swift, and H. L. Swinney, Phys. Rev. 
Lett. 75, 4397 (1995). 

6 U. Thiclc, M. Mcrtig, and W. Pompe, Phys. Rev. Lett. 80, 2869 (1998). 

7 A. Oron, S. H. Davis, and S. G. Bankoff, Rev. Mod. Phys. 69, 931 (1997). 
^ O. Reynolds, Phil. Trans. Roy. Soc. 177, 157 (1886). 

9 A. Sommerfeld, Z. Math. Phys. 50, 97 (1904). 

° A. Sharma and G. Reiter, J. Colloid Interface Sci. 178, 383 (1996). 

^ G. Reiter, Langmuir 9, 1344 (1993). 

^ K. Jacobs, R. Seemann, G. Schatz, and S. Herminghaus, Langmuir 14, 4961 (1998). 
3 S. J. VanHook, M. F. Schatz, J. B. Swift, W. D. McCormick, and H. L. Swinney, J. Fluid Mech. 
345, 45 (1997). 

^ R. Seemann, S. Herminghaus, and K. Jacobs, Phys. Rev. Lett. 86, 5534 (2001). 

^ E. Ruckenstein and R. K. Jain, J. Chcni. Soc. Faraday Trans. II 70, 132 (1974). 

6 M. B. Wilhams and S. H. Davis, J. Colloid Interface Sci. 90, 220 (1982). 

^ P. G. de Gennes, Rev. Mod. Phys. 57, 827 (1985). 

^ S. G. Yiantsios and B. G. Higgins, Phys. Fluids A 1, 1484 (1989). 

^ A. Oron and P. Rosenau, J. Physique II France 2, 131 (1992). 

'0 A. Sharma, Langmuir 9, 861 (1993). 

A. Sharma, Langmuir 9, 3580 (1993). 
'2 V. S. Mitlin, J. Colloid Interface Sci. 156, 491 (1993). 

A. Sharma and R. Khanna, Phys. Rev. Lett. 81, 3463 (1998). 
'4 U. Thiele, M. G. Velarde, and K. Neuffer, Phys. Rev. Lett. 87, 016104 (2001). 

E. Schaffer and U. Steiner, Eur. Phys. J. E 8, 347 (2002). 
'6 U. Thiele, K. Neuffer, Y. Pomeau, and M. G. Velarde, Colloid Surf. A 206, 135 (2002). 

U. Thiele, Eur. Phys. J. E 12, 409 (2003). 



24 



2^ A. L. Bertozzi, G. Griin, and T. P. Witelski, Nonlinearity 14, 1569 (2001). 

29 M. Bestehorn, A. Pototsky, and U. Thiele, Eur. Phys. J. B 33, 457 (2003). 

3° U. Thiele and E. Knobloch, Physica D 190, 213 (2004). 

^1 A. Pototsky, M. Bestehorn, and U. Thiele, Physica D 199, 138 (2004). 

^2 A. A. Golovin, A. A. Nepomnyashchy, and L. M. Pismen, Phys. Fluids 6, 34 (1994). 

J. N. Israelachvili, Inter-molecular and Surface Forces (Academic Press, London, 1992). 

R. J. Hunter, Foundation of Colloid Science, vol. 1 (Clarendon Press, Oxford, 1992). 
^5 R. Seemann, S. Herminghaus, and K. Jacobs, J. Phys.-Cond. Mat. 13, 4925 (2001). 
36 Z. Q. Lin, T. Kerle, S. M. Baker, D. A. Hoagland, E. Schaffer, U. Steiner, and T. P. Russeh, J. 

Chem. Phys. 114, 2377 (2001). 
3^ Z. Q. Lin, T. Kerle, T. P. Russell, E. Schaffer, and U. Steiner, Macromolecules 35, 3971 (2002). 

38 D. Merkt, A. Pototsky, M. Bestehorn, and U. Thiele, Phys. Fluids (2005), (submitted). 

39 F. Brochard-Wyart, P. Martin, and C. Redon, Langmuir 9, 3682 (1993). 

'^^ Y. L. Zhang, O. K. Matar, and R. V. Craster, J. Colloid Interface Sci. 262, 130 (2003). 
^1 R. V. Craster and O. K. Matar, J. Fluid Mech. 425, 235 (2000). 
^2 O. K. Matar, R. V. Craster, and M. R. E. Warner, J. Fluid Mech. 466, 85 (2002). 
^3 K. D. Danov, V. N. Paunov, N. Alleborn, H. Raszillier, and F. Durst, Chem. Eng. Sci. 53, 2809 
(1998). 

K. D. Danov, V. N. Paunov, S. D. Stoyanov, N. Alleborn, H. Raszillier, and F. Durst, Chem. 
Eng. Sci. 53, 2823 (1998). 

V. N. Paunov, K. D. Danov, N. Alleborn, H. Raszillier, and F. Durst, Chem. Eng. Sci. 53, 2839 
(1998). 

^6 D. Bandyopadhyay, Stability and dynamics of bilayers (2001), master-thesis. Dep. Chem. Eng., 

Ind. Inst. Tech. Kanpur. 
^7 A. Pototsky, M. Bestehorn, D. Merkt, and U. Thiele, Phys. Rev. E 70, 025201 (R) (2004). 
^8 O. Wunnicke, P. Miiller-Buschbaum, M. Wolkenhauer, C. Lorenz-Haas, R. Cubitt, V. Leiner, 

and M. Stamm, Langmuir 19, 8511 (2003). 
'^^ M. D. Morariu, E. Schaffer, and U. Steiner, Phys. Rev. Lett. 92, 156102 (2004). 
^° A. Faldi, R. J. Composto, and K. I. Winey, Langmuir 11, 4855 (1995). 
5^ P. Lambooy, K. C. Phelan, O. Haugg, and G. Krausch, Phys. Rev. Lett. 76, 1110 (1996). 
52 Q. Pan, K. I. Winey, H. H. Hu, and R. J. Composto, Langmuir 13, 1758 (1997). 

25 



C. Renger, P. Muller-Buschbaum, M. Stamm, and G. Hinrichsen, Macromolecules 33, 8388 
(2000). 

5^ P. Martin, A. Buguin, and F. Brochard-Wyart, Europhys. Lett. 28, 421 (1994). 

A. A. Nepomnyashchii and I. B. Simanovskii, Pmm-J. Appl. Math. Mech. 54, 490 (1990). 

A. A. Nepomnyashchy and I. B. Simanovskii, Q. J. Mech. Appl. Math. 50, 149 (1997). 
5^ G. F. Teletzke, H. T. Davis, and L. E. Scriven, Rev. Phys. Appl. 23, 989 (1988). 

B. V. Deryagin and L. D. Landau, Acta Physicochim. 14, 633 (1941). 

E. J. W. Verwey and J. T. G. Overbeck, Theory of the Stability of Lyophillic Colloids (Elsevier, 
Amsterdam, 1948). 

60 H. Ohshima, Colloid Polym. Sci. 252, 158 (1974). 

61 H. Ohshima, Colloid Polym. Sci. 252, 257 (1974). 

6^ E. S. Manteuffel and K. Vetters, Lineare Algebra (Teubner Verlagsgesellschaft, Leipzig, 1972). 

63 M. D. Morariu, E. Schaffer, and U. Steiner, Eur. Phys. J. E 12, 375 (2003). 

64 M. O. David, G. Reiter, T. Sitthai, and J. Schultz, Langmuir 14, 5667 (1998). 

65 M. Sferrazza, C. Xiao, R. A. L. Jones, D. G. Bucknall, J. Webster, and J. Penfold, Phys. Rev. 
Lett. 78, 3693 (1997). 

66 S. G. Yiantsios and B. G. Higgins, J. Colloid Interface Sci. 147, 341 (1991). 
6^ E. Docdcl, H. B. Keller, and J. P. Kerncvez, Int. J. Bif. Chaos 1, 493 (1991). 
6^ E. Doedel, H. B. Keller, and J. P. Kernevez, Int. J. Bif. Chaos 1, 745 (1991). 

6^ E. Doedel, A. Champneys, T. Fairfrieve, Y. Kusnetzov, B. Sandstede, and X. Wang, Con- 
tinuation and Bifurcation Software for Ordinary Differential Equations (Concordia University, 
Montreal, 1997). 

^0 U. Thiele, L. Brusch, M. Bestehorn, and M. Bar, Eur. Phys. J. E 11, 255 (2003). 
L. Mahadevan, M. Adda-Bedia, and Y. Pomeau, J. Fluid Mech. 451, 411 (2002). 



26 



TABLE I: Possible combinations of signs of the different Hamaker constants for given order of the 
refractive indices of the involved medias. 



refractive indices 


^129 


A21S 




Ag21s 


Us > ni, ni < 712, n2 > rig 


+ 


+ 






Us < ni, ni < n2, n2 > Ug 


+ 






+ 


Hg < ni, ni > n2, n2 > Ug 




+ 




+ 


Ug > 111, rii > 712, n2> Ug 


TABLE IL Hamaker constants for various combinations of polymers. 


System 


Ai2g X 10-20Nm 


A21S X 10-20Nm 


Ag21s 


X 10-20Nm 


Si/PMMA/PS/air 


1.49 


3.8 




-23.02 


SiO/PMM A/PS/air 


1.49 


-0.024 




0.15 


SiO/PS/PDMS/air 


-1.83 


0.42 




1.25 
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FIG. 1: Sketch of the problem in two dimensions. The local thickness of the lower layer is hi, the 
total local film thickness is h2- 

FIG. 2: The schematic stability diagram for fixed coupling E12. Shown are the stability threshold 
(solid line) and the boundary between the one- mode and the two-mode region (dotted line). Both 
are given by det E = for increasing wave number k. Dashed and dot-dashed arrows represent 
parametric lines given by {Eii{k) and E22{k)). The dashed (dot-dashed) line starts at /c = in 
the two-mode (one-mode) region. At an intersection of a line det E{k = 0) and a parametric line 
one of the growth rates changes its sign. 

FIG. 3: Shown are the growth rate 7 (solid lines) and the mode type (j) (dashed lines) of the leading 
eigenmode. (a) A varicose mode from the one-mode region at di = 30, d2 = 47 and a = fi = 1, 
(b) a zigzag mode from the one- mode region at di = 15, d2 = 40 and a = fi = 1. Panel (c) gives 7 
and (p for di and d2 as in (b) but for /x = 0.1. For convenience we plot in (b) IO7 and in (a) 2O7. 

FIG. 4: Different types of stability diagrams in the plane of the layer thicknesses {di,d2 — di), 
shown for different strength of the short-range interactions Si and S2 as given in the legends. 
The shaded parts represent linearly stable regions. The Hamaker constants are ^i2g = 1-49, 
A21S = 3.8, Ag2is = -23.02, corresponding to the Si/PMMA/PS/air system. Panels (a), (b), (c) 
and (d) correspond to ranges I, III, IV and II in Fig. [51 respectively. 

FIG. 5: Phase diagram in the plane (5i, 5*2) for the Si/PMMA/PS/air system. The absolute 
stable region (hatched rectangle in the lower left corner) is bounded by (5*1)111111 = (e/4)^^2is/|^i2g| 
from the right and by (5'2)min = (e/4)^ from above. The unstable region is divided into seven 
qualitatively different subregions, described in the main text. 

FIG. 6: Snapshots of the time evolution of a Si/PMMA/PS/air system for di = 30, d2 = 47, 
5*1 = 5*2 = 1, cj = 1 and /x = 1 at times as given in the legends. The domain length is L = AX^ and 
time is in units of l/7m- 
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FIG. 7: Evolution in time of (a) the relative energy E — Eq and (b) the integral mode-type (pint 
fEa. l32() for parameters as in Fig. |H1 Time is in units of l/'jm- In (a) the dashed lines denote the 
energy levels which correspond to the stationary solutions with periods L = Am (first and second 
line from above), L = 4/3Am (third line from above) and L = 2 Am (the lowest line), taken from 
Fig.Hb). 

FIG. 8: Characterization of the stationary periodic solutions for the system of Fig.[7| Shown are 
(a) the amplitude of the lower layer Ai, (b) the relative energy E — Eq, (c) the norm L2, and (d) 
the integral mode type (pint in their dependence on the period L. In (d) the inset shows a zoom of 
the region marked by the dashed box. 

FIG. 9: The two stationary solutions with period L = 108.28 (cp. Fig.|HI) corresponding to the 
wave length of the relevant dominant linear mode Xmax- Shown are (a) the varicose type and (b) 
the zigzag type solution from the first and second stable branch, respectively. 

FIG. 10: Snapshots of the time evolution of a Si/PMMA/PS/air system for di = 15, d2 = 40, 
5*1 = 5*2 = 1, o" = 1 and = 1 at times as given in the legends. The domain length is L = 4 Am 
and time is in units of l/7m- 

FIG. 11: Evolution in time of (a) the relative energy E — Eq and (b) the integral mode- type (pint 
(Eq. 132(1 for parameters as in Fig. 1101 The inset in (b) shows the early-time behaviour of the mode- 
type. Time is in units of l/7m- In (a) the dashed lines denote the energy levels which correspond 
to the stationary solutions with periods L = Am, L = 4/3Am and L = 2 Am, taken from Fig. ll2r b'l. 

FIG. 12: Characterization of the stationary periodic solutions, for the system of Fig. llUI Shown 
are (a) the amplitude of the upper layer A2, (b) the relative energy E — Eq, (c) the norm L2 and 
(d) the integral mode type (pint in their dependence on the period L. 

FIG. 13: The two stationary solutions with period L = 50.8 and L = 101.6 (cp. Fig. 112(1 corre- 
sponding to once and twice the wave length of the dominant linear mode Amax, respectively. To 
symbolize the coarsening process we show in (a) and (b) two and one period(s), respectively. The 
X coordinate is in units of Xmax- 
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FIG. 14: Large-period (long-time) stationary profiles for (a) a Si/PMMA/PS/air system with 
di = 30, d2 = 39, Xm = 132 and period L = 26 x A^, (b) a SiO/PMMA/PS/air with di = 30, 
d2 = 50, = 246 and period L = 20 x A„, (c) a Si/PMMA/PS/air system with di = 30, d2 = 70, 
Xm = 118 and L = 21 x A^ and (d) a SiO/PS/PDMS/air with di = 30, d2 = 70, A^ = 557 and 
L = 10 X A^. The remaining parameters are cr = 1, ;u = 1, and Si = S2 = 1- 



FIG. 15: Sketch of the system with a gas layer of finite thickness. 
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